Analysis of stimulation experiments.
library(Seurat)
## Attaching SeuratObject
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(ggplot2)
Load integrated object. This object has been integrated with k.anchor = 20, because there was lots of batch effect between the 0hr tissue and the other two experiments (which merged quite well).
load("~/Dropbox/Zoe/scf_version/analysis/healthy_sc/seurat_objects/no_dropletQC/3391/integrated_3391_stimulation_kanchor20_scClustViz.RData")
Set identities:
Idents(scSeurat) <- "integrated_snn_res.0.4"
Check out the UMAP of clusters, original identities, and refined SCINA labels.
DimPlot(scSeurat, label = TRUE)
DimPlot(scSeurat, group.by = "tissue_type", label = FALSE)
DimPlot(scSeurat, group.by = "scina_labels_refined", label = TRUE) & NoLegend()
Look at tissue-type per cluster on a proportional level:
pt <- table(Idents(scSeurat), scSeurat$tissue_type)
pt <- as.data.frame(pt)
ggplot(pt, aes(x = Var1, y = Freq, fill = Var2)) +
theme_bw(base_size = 15) +
geom_col(position = "fill", width = 0.5) +
xlab("Sample") +
ylab("Proportion") +
theme(legend.title = element_blank())
Look for T-cell signatures.
FeaturePlot(scSeurat, features = c("CD3E", "NKG7"))
FeaturePlot(scSeurat, features = c("LEF1", "XCL1;XCL2"))
FeaturePlot(scSeurat, features = c("TOX", "CD8A"))
Look for macrophage signatures:
FeaturePlot(scSeurat, features = c("LGALS1", "TYROBP"))
FeaturePlot(scSeurat, features = c("MARCO", "C1QC"))
FeaturePlot(scSeurat, features = c("LYZ-1", "S100A8"))
FeaturePlot(scSeurat, features = c("DEFA6;DEFA4;DEFA5", "CTSS"))
FeaturePlot(scSeurat, features = c("FLT3", "S100A12"))
Cluster 11 is the only cluster that is obviously dominated by the 0hr and unstimulated experiments. What cells are these? Looks like they are hepatocytes. A lot of the markers are periportal, but all hepatocytes are in this cluster.
# Isolate DE genes
clust11 <- sCVdata_list$res.0.4@DEvsRest$`11`[order(sCVdata_list$res.0.4@DEvsRest$`11`$Wstat,
decreasing = TRUE), ]
rownames(clust11)[1:50]
## [1] "PAH"
## [2] "LOC114082496"
## [3] "GYS2"
## [4] "SLC2A2"
## [5] "SLCO1B3;SLCO1B3-SLCO1B7;SLCO1B7;SLCO1B1"
## [6] "PCK1"
## [7] "CPS1"
## [8] "AGMO"
## [9] "MTSS1"
## [10] "FBP1"
## [11] "ALDOB"
## [12] "MAT1A"
## [13] "mikado.WCK01-AAH20201022-F8-rsc00143G668"
## [14] "SHMT1"
## [15] "CYP3A4;CYP3A7-CYP3A51P;CYP3A7-3"
## [16] "KLKB1"
## [17] "TDO2"
## [18] "SLC38A4"
## [19] "SORBS2"
## [20] "SOX5"
## [21] "CYP3A4;CYP3A7-CYP3A51P;CYP3A7-2"
## [22] "ARG1"
## [23] "SLC25A13"
## [24] "CYP3A4;CYP3A7-CYP3A51P;CYP3A7"
## [25] "Slc22a19-1"
## [26] "TAT"
## [27] "LOC107139912"
## [28] "KIRREL3"
## [29] "ACACB"
## [30] "C5-1"
## [31] "SUGCT"
## [32] "ZBTB16"
## [33] "ABCB11"
## [34] "CFH-1"
## [35] "GAS2"
## [36] "ACSL5"
## [37] "STARD13"
## [38] "IGFBP1"
## [39] "HSD11B1"
## [40] "AQP9"
## [41] "GRHPR"
## [42] "FMO5-1"
## [43] "CMBL"
## [44] "CYP7B1"
## [45] "C8G"
## [46] "IDO2"
## [47] "SLC41A2"
## [48] "DMD"
## [49] "mikado.WCK01-AAH20201022-F8-rsc00002G194"
## [50] "LOC107139914"
FeaturePlot(scSeurat, features = c("PAH","PCK1"))
FeaturePlot(scSeurat, features = c("CYP2E1","FETUB"))
FeaturePlot(scSeurat, features = c("ACACB","ELOVL6"))
FeaturePlot(scSeurat, features = c("POLR2D", "ALB-1"))
Let’s look at some of the key markers that Sonya sent: IL2, IL4, IL6, IL7, IL10, IL12A, TGFB
FeaturePlot(scSeurat, features = c("IL2", "IL4"))
FeaturePlot(scSeurat, features = c("IL6", "IL6-1"))
FeaturePlot(scSeurat, features = c("IL7", "IL10"))
FeaturePlot(scSeurat, features = c("IL12A", "TGFB1"))
It’s pretty hard to see which cytokines are expressed in which samples, so I’m going to try splitting up the samples.
DimPlot(scSeurat, split.by = "tissue_type", label = TRUE) & NoLegend()
Now let’s take a look at the features on this split plot.
FeaturePlot(scSeurat, feature = "IL2", split.by = "tissue_type") # Only in PMAIO cluster 1
## Warning in FeaturePlot(scSeurat, feature = "IL2", split.by = "tissue_type"):
## All cells have the same value (0) of IL2.
## Warning in FeaturePlot(scSeurat, feature = "IL2", split.by = "tissue_type"):
## All cells have the same value (0) of IL2.
FeaturePlot(scSeurat, feature = "IL4", split.by = "tissue_type") # Only in PMAIO cluster 1
## Warning in FeaturePlot(scSeurat, feature = "IL4", split.by = "tissue_type"):
## All cells have the same value (0) of IL4.
## Warning in FeaturePlot(scSeurat, feature = "IL4", split.by = "tissue_type"):
## All cells have the same value (0) of IL4.
FeaturePlot(scSeurat, feature = "IL6", split.by = "tissue_type") # Stronger in endo and mesem PMAIO
## Warning in FeaturePlot(scSeurat, feature = "IL6", split.by = "tissue_type"):
## All cells have the same value (0) of IL6.
FeaturePlot(scSeurat, feature = "IL6-1", split.by = "tissue_type") # No signal
## Warning in FeaturePlot(scSeurat, feature = "IL6-1", split.by = "tissue_type"):
## All cells have the same value (0) of IL6-1.
## Warning in FeaturePlot(scSeurat, feature = "IL6-1", split.by = "tissue_type"):
## All cells have the same value (0) of IL6-1.
FeaturePlot(scSeurat, feature = "IL7", split.by = "tissue_type") # More in PMAIO cluster 4
FeaturePlot(scSeurat, feature = "IL10", split.by = "tissue_type") # More in PMAIO cluster 1
## Warning in FeaturePlot(scSeurat, feature = "IL10", split.by = "tissue_type"):
## All cells have the same value (0) of IL10.
FeaturePlot(scSeurat, feature = "IL12A", split.by = "tissue_type") # No signal
## Warning in FeaturePlot(scSeurat, feature = "IL12A", split.by = "tissue_type"):
## All cells have the same value (0) of IL12A.
FeaturePlot(scSeurat, feature = "TGFB1", split.by = "tissue_type") # Pretty even across all
Cluster 1 looks like a really interesting cluster that has really expanded in the PMAIO experiment. Isolate DE genes for this cluster.
clust1 <- sCVdata_list$res.0.4@DEvsRest$`11`[order(sCVdata_list$res.0.4@DEvsRest$`1`$Wstat,
decreasing = TRUE), ]
cat(rownames(clust1)[1:50])
## MIB1 FBXL4 PPP3CA OPA1 ST7L CYRIB FILIP1L ARG2 RARB P3H1 TCF4 VPS41 OGA SLC19A2 LOC114104752 AFDN SRP9 ABCC2 TBC1D19 PANK1 FBXO34 ARIH2 OSBPL9-1 SLC22A3 AZGP1 MSH3 SPP2 FBXO25 DLEC1 FPGS TLE4 SLC10A7 B2M HDAC1 TRAM1 FAM135A ANGPTL3 MKKS ABLIM3 TTLL5 B4GALT1 TBC1D12 FKBP8 SCLY PDHX EIF4G2 CNOT10 mikado.WCK01-AAH20201022-F8-rsc00035G807 PACS2 CEP350
Visualise these genes in violin plots.
Interpretation: Now I see that IL2 and IL4 are expressed in more than just cluster 1, but it is nearly always in PMAIO. Cluster 8 is definitely where IL6 is expressed. IL7 is a bit more spread out across treatment types, but it’s definitely present in PMAIO in clusters 4, 5, and 14. Cluster 14 is macrophages, but I’m honestly not sure of 4 and 5 - low signal hepatocytes? IL10 is in clusters 1 and 13 PMAIO. IL12A is a bit more in 6 PMAIO.
VlnPlot(scSeurat, feature = "IL2", split.by = "tissue_type") # Only in PMAIO cluster 1
## The default behaviour of split.by has changed.
## Separate violin plots are now plotted side-by-side.
## To restore the old behaviour of a single split violin,
## set split.plot = TRUE.
##
## This message will be shown once per session.
VlnPlot(scSeurat, feature = "IL4", split.by = "tissue_type") # Only in PMAIO cluster 1
VlnPlot(scSeurat, feature = "IL6", split.by = "tissue_type") # Stronger in endo and mesem PMAIO
VlnPlot(scSeurat, feature = "IL6-1", split.by = "tissue_type") # No signal
VlnPlot(scSeurat, feature = "IL7", split.by = "tissue_type") # More in PMAIO cluster 4
VlnPlot(scSeurat, feature = "IL10", split.by = "tissue_type") # More in PMAIO cluster 1
VlnPlot(scSeurat, feature = "IL12A", split.by = "tissue_type") # No signal
VlnPlot(scSeurat, feature = "TGFB1", split.by = "tissue_type") # Pretty even across all